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We propose a general method for the construction and analysis of unweighted t - recurrence 
networks from chaotic time series. The selection of the critical threshold e c in our scheme is done 
empirically and we show that its value is closely linked to the embedding dimension M. In fact, 
we are able to identify a small critical range Ae numerically that is approximately the same for the 
random and several standard chaotic time series for a fixed M. This provides us a uniform framework 
for the non subjective comparison of the statistical measures of the recurrence networks constructed 
from various chaotic attractors. We explicitly show that the degree distribution of the recurrence 
network constructed by our scheme is characteristic to the structure of the attractor and display 
statistical scale invariance with respect to increase in the number of nodes N. We also present two 
practical applications of the scheme, detection of transition between two dynamical regimes in a time 
delayed system and identification of the dimensionality of the underlying system from real world 
data with limited number of points, through recurrence network measures. The merits, limitations 
and the potential applications of the proposed method have also been highlighted. 

PACS numbers: 05.45.Ac, 05.45.Tp, 05.45.Df 


I. INTRODUCTION 

Study of networks has become an important area of 
research in the last two decades cutting across various 
disciplines and often providing a coherent view of struc¬ 
tures and phenomena which may appear different from a 
common perspective. Mathematically, networks are en¬ 
tities defined on an abstract space, with N number of 
nodes and arbitrary number of links between them. We 
often refer to complex networks , implying that the struc¬ 
ture is irregular and complex. Such networks may be 
weighted or unweighted, directed or undirected depend¬ 
ing on the structure and interaction of the system it tries 
to model. 

Historically, the study of networks has been the do¬ 
main of a branch of mathematics called graph theory. 
The mathematical basis for the analysis of complex net¬ 
works was laid using the so called “random graphs” (RG) 
by Erdos and Renyi about five decades back |1|. They 
have shown several properties for random graphs, which 
are now popularly known as Erdos-Renyi (E-R) networks. 
The measures introduced for the analysis of RGs have, 
naturally, become the tools to characterize the structure 
and evolution of complex networks. The most impor¬ 
tant among them are the degree distribution and the 
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characteristic path length (GPL), apart from other mea¬ 
sures, such as, link density (LD), clustering coefficient 
(CC), diameter, centrality, etc. The degree distribution 
indicates how many nodes nk among the total number 
of nodes N have a given degree k. It is usually repre¬ 
sented as a probability distribution P{k ) as a function of 
k, where P{k) = . For E-R networks, one can show 

that P(k) follows a Binomial distribution which, for large 
N and small probability of connection tends to the Pois¬ 
son distribution. The CPL, denoted by < l >, is defined 
through the shortest path l s connecting two nodes i and 
j. For unweighted and undirected networks that we con¬ 
sider in this work, l s is defined as the minimum number 
of nodes to be traversed to reach from i to j. The average 
value of l s for all the pair of nodes in the whole network 
is defined as < l > and the maximum value of l s is taken 
as the diameter of the network, denoted by U. For a 
detailed discussion of all the network measures, see the 
popular books by Newman [2] and Watts fU and some 
excellent reviews on the subject Hi- 

An important area where the new network based con¬ 
cepts and measures have been applied successfully is in 
the analysis of dynamical systems @, [7], especially those 
showing chaotic behavior, by constructing complex net¬ 
works from time series data of the dynamical systems. 
Here, the aim is to extract information regarding the 
structure of the underlying chaotic attractor, which are 
otherwise difficult to get using the conventional methods 
of nonlinear time series analysis. The basic idea of this 
technique is that the information inherent in a chaotic 
time series is mapped on to the domain of a complex 
network using a suitable scheme. One then uses the sta- 
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tistical measures of the complex network to characterize 
the underlying chaotic attractor. 

Two questions are relevant in this context. Firstly, 
which method is to be used to transform the time se¬ 
ries into the corresponding network and secondly, how 
to ensure that the resulting network truly represents the 
characteristic features of the underlying attractor. To 
answer the first, several methods have been suggested 
in the literature, such as cycle networks {s], visibility 
graphs transition networks 0 and recurrence net¬ 
works (RN) & The RNs can be constructed in differ¬ 
ent ways, namely, the correlation networks HU . ^-nearest 
neighbour networks ITdf and e-recurrence networks 141. 
It has been shown that the networks generated by each of 
these methods can capture several characteristics of the 
chaotic time series like dynamical transitions in the sys¬ 
tem, topological properties of the attractor, etc. Hence 
each method is relevant in the context of specific appli¬ 
cations. A detailed discussion of these methods and their 
comparison can be found elsewhere 01- 

Among the methods mentioned above, the one based 
on e - recurrence is physically appealing since it is based 
on the concept of recurrence of a trajectory in the phase 
space. Moreover, the method can be applied to any type 
of synthetic or real world data and the resulting networks 
are found to be useful tools for uncovering complex bi¬ 
furcation scenario and detecting dynamical transitions 
in palaeo-climate data EMI- The methods based on 
recurrence have also found diverse applications ranging 
from life sciences 0,0], earth sciences 221 and astro¬ 
physics [23j. In this work, we concentrate on this method 
for network generation, confining ourselves to the case of 
unweighted e-recurrence networks. 

The answer to the second question raised above leads 
us to the choice of the parameters for network construc¬ 
tion. A crucial parameter in the construction of the RN 
is the recurrence threshold e itself. In all the existing 
schemes, the value of e chosen is different for each time 
series, whatever be the criteria used for its selection, due 
to the arbitrary size of the attractor after embedding. 
Our main goal in the present work is to propose a scheme 
that uses an approximately identical range of values for 
e for different time series, both synthetic and real world, 
for a given embedding dimension. Apart from providing 
a uniform framework for the recurrence network analysis, 
the scheme has several advantages and practical applica¬ 
tions as discussed below in detail. 


Our paper is organised as follows: In the next section, 
we discuss the basic idea of RN construction while in §111, 
the criteria for the selection of all the parameters for the 
construction of RN from the time series are presented. 
We then proceed, in §IV, to construct the RNs from sev¬ 
eral low dimensional chaotic systems. All the important 
network measures are derived from the RNs as a function 
of M and N and compared. The degree distribution, es¬ 
pecially, is studied in detail and is shown to be character¬ 
istic of the structural complexity of a chaotic attractor. 
Two practical applications of the proposed scheme are 


illustrated in §V. A discussion on various aspects of im¬ 
plementation of the scheme and conclusions are given in 
§VI. 


II. CONSTRUCTION OF RECURRENCE 
NETWORK 

In this section, we briefly discuss the basic idea regard¬ 
ing the construction of recurrence networks. More de¬ 
tails can be found in recent reviews on the topic flU ITfl . 
Recurrence is a fundamental property of every bounded 
dynamical system by which a trajectory tends to revisit 
a certain region of the phase space over a time inter¬ 
val. This basic concept has been utilized to develop a 
visualization tool called the recurrence plot (RP) for the 
analysis of dynamical systems [24|. A RP represents all 
recurrences in the form of a binary matrix 1Z where Rij 
= 1 if the state Xj is a neighbour of Xi in phase space 
and Rij = 0, otherwise. The neighbourhood is defined 
through a certain recurrence threshold e. In the most 
general definition, the discretely sampled scalar time se¬ 
ries s(l),s(2),. s(Nt) is embedded in M-dimensional 

space taking the time delay co-ordinates [2b[ using a suit¬ 
able time delay r, where Nt is the total number of points 
in the time series. The procedure creates delay vectors 
in the embedded space of dimension M given by 

Xi = [s(«), s(i + t),. s(i + (M - l)r)] (1) 

There are a total number of N = Nt—(M—1)t state vec¬ 
tors in the reconstructed space representing the attrac¬ 
tor. Any point j on the attractor is considered to be in 
the neighbourhood of a reference point 1 if their distance 
in the M-dimensional space is less than the threshold e. 
Thus we have 


Rij = H(e-\\xi-Xj\\) (2) 


where H is the Heaviside function and ||..|| is a suitable 
norm. In this paper, we use the Euclidean norm. The RP 
can only visually distinguish between different qualitative 
features of dynamics. This tool has become more popu¬ 
lar with the introduction of the recurrence quantification 
analysis (RQA) 0] using the measures derived from the 
RP. It has found numerous applications 0-0 and even 
dynamical invariants like correlation dimension D 2 and 
correlation entropy I \2 can be evaluated efficiently using 
RQA 0. 

The importance of the e - RN (which, from now on, we 
simply call RN) is that its generation is closely associated 
with the RP. In fact, the adjacency matrix A for the 
unweighted RN can be obtained by removing the identity 
matrix from the recurrence matrix: 

A'ij — Rij 3ij ( 3 ) 

where Sij is the Kronecker delta. Note that, once the 
adjacency matrix is defined, the time series has been 
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converted into a complex network. Each point on the 
embedded attractor is taken as a node in the RN and a 
node i is connected to another node j if the distance dij 
between the corresponding points on the embedded at¬ 
tractor is < e. Thus the adjacency matrix A is a binary 
symmetric matrix with elements Aj 3 = 1 if < e and 
0 otherwise. Note that, in contrast to the RP measures 
which consider the temporal properties of the trajectory 
points, RN analysis quantifies the geometrical properties 
of the underlying attractor and hence can give useful in¬ 
formation regarding the structure of the attractor. 

Even though the method for the generation of the RN 
appears to be simple, it has several ambiguities associ¬ 
ated with it (33J. How do we ensure that the RN captures 
the structural characteristics of the underlying attractor? 
The answer lies in the proper choice of the parameters 
involved in the construction of the network. For the RN, 
the key parameters are e and M. If e is large, specific 
small scale properties of the attractor cannot be revealed 
and if e is too small, the network breaks into dissuaded 
nodes due to lack of connections. Many authors dlEl- 
[3l have discussed this issue in detail and have given some 
guidelines for the choice of e. But for arbitrary size of the 
attractor, the choice of e still remains subjective. Simi¬ 
larly, the specific feature of RN generation is embedding 
and in this context the choice of M has not been dis¬ 
cussed much in the literature as it is commonly believed 
that M should be sufficiently high for the attractor to be 
fully resolved. We show that the choice of e is closely re¬ 
lated to that of M and we present a scheme for the choice 
of e and M that gives a uniform critical range of e for all 
time series for a given M. To validate the wide range 
of applicability of the scheme, we show results from sev¬ 
eral low dimensional chaotic attractors as well as random 
data. 


III. CHOICE OF PARAMETERS FOR 
NETWORK CONSTRUCTION 

There are four parameters associated with the RN gen¬ 
eration, which are the time delay r, e, M and N, the 
number of nodes. Note that N < Nt , the total number 
of points in the time series and the difference depends 
on M and r. The value of Nt can be adjusted to get 
the required number of N for the computation. For the 
choice of r we stick to the most commonly used criteria, 
namely, the first minimum of the autocorrelation func¬ 
tion. The value of r is related to the time step At used 
for the generation of the time series. For the sake of uni¬ 
formity, we use At = 0.05 to generate the time series 
from all the continuous time systems presented here. We 
have removed the first 10000 values as transients in all 
cases. 

In all our numerical computations we use the value of 
N in the range 2000 to 10000. The lower limit is set be¬ 
cause, if the number of data points in the time series is 
too small, the basic structure of the embedded attractor- 


may not have evolved completely. The upper limit is set 
mainly due to the fact that the computations become 
increasingly difficult due to high memory requirement 
for N > 10000. However, we find that N < 10000 is 
sufficient to get reasonable results from low dimensional 
chaotic systems. Moreover, for many real world applica¬ 
tions of RN analysis, one has to confine to this range of 
N very often. 

We next consider the choice of the crucial parameter, 
namely, the critical threshold, e c . There are already a 
number of criteria suggested for choosing e c . For ex¬ 
ample, Gao and Jin [321 ] gives a heuristic criterion that 
selects e c as the value at which the link density becomes 
maximum when plotted as a function of e. But this 
has the drawback that small changes in e induce large 
changes in the network measures as indicated by Donner 
et al. [3l]]. Recently, another criterion has been sug¬ 
gested based on analytic methods [34] while an adaptive 
selection of threshold has been proposed by Eroglu et al. 
for specific time series [35j. The last one is especially 
important as it chooses the critical threshold based on 
the theory of random graphs, where the second small¬ 
est eigen value A 2 of the Laplacian matrix C crosses zero 
when plotted as a function of e as the network becomes 
fully connected [sj, where C = V — A, the difference be¬ 
tween diagonal degree matrix and adjacency matrix. We 
will show below that this criterion comes very close to 
the empirical criterion used by us. 

In all the above works so far considered, the size of the 
attractor after embedding is arbitrary so that the value 
of the threshold will be different for different attractors. 
Our primary motivation in the present work is to look for 
a scheme that can give approximately identical value for 
the critical threshold for different time series. We con¬ 
sider this to be an important step forward as it will lead 
to a uniform framework for the RN analysis which may 
be more useful for application to practical time series, as 
shown below. 

To overcome the problem of arbitrary size of the attrac¬ 
tor, we transform the time series into a uniform deviate. 
For this, we first rescale the time series into the unit in¬ 
terval [0,1]. We then take each value yi in the time series 
and count how many values are less than or equal to yi. 
Let this count be rq. Then the uniform deviate time 
series Ui is obtained as Ui = where N T is the total 
number of points in the time series. The effect of uniform 
deviate transformation is shown in Fig. [l] for the Lorenz 
attractor, where the original time series y(t) and the time 
series after uniform deviate u(t) are shown along with 
the corresponding attractors after embedding. We have 
shown the importance of uniform deviate transformation 
in computing the conventional nonlinear measures like 
correlation dimension D 2 and entropy K 2 HM , espe- 
daily from higher dimensional systems [Hj]. It stretches 
the embedded attractor uniformly in all directions with¬ 
out changing any of the dynamical invariants of the at¬ 
tractor, minimises the edge effect and provides improved 
scaling region and better convergence with data points. 
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FIG. 1: (Color online) Top panel shows the original time 
series y(t) from the Lorenz attractor on the left and the time 
series u(t) after uniform deviate transformation on the right. 
The bottom panel shows the corresponding attractors after 
embedding. 




The primary criterion that we use for the selection of 
e c in our scheme is that the resulting RN has to remain 
mostly as “one single cluster”. This is checked empiri¬ 
cally by changing e. Note that this criterion is analogous 
to the one suggested by Donges et al. [HJ where the au¬ 
thors propose to choose the percolation threshold above 
which the giant component of the RN appears. For the 
selection of e c , we use a more practical approach rather 
than the rigorous criterion as adopted by the previously 
mentioned authors and we show the advantages of our 
approach in the sections below. A formal method to se¬ 
lect e c is to choose the value of e when the network just 
becomes fully connected. This is done by computing the 
second smallest eigen value A 2 of the Laplacian matrix 
as a function of e. From the theory of random graphs, it 
is well known jhj that as A 2 crosses zero from negative, 
the network becomes fully connected, with no dissuaded 
nodes. In Fig. [2J we show the variation of A 2 with e for 
RNs constructed from Lorenz and random time series, 
for N = 2000 and 5000, with M = 3. It is evident that 
A 2 for random network becomes positive slightly earlier 
compared to Lorenz in both cases. As N increases, there 
is also a slight shift towards lower e. The network for 
random becomes fully connected for e = 0.09 while for 
Lorenz, the value is 0.13. We have repeated the computa¬ 
tion for other standard chaotic attractors and found that 
the value of e where A 2 becomes positive varies slightly 
for different systems. However, we find that there is a 
small range of e, from 0.1 to 0.13, where the network is 
almost fully connected for all systems with the appear¬ 
ance of a giant cluster and very few (< 1%) dissuaded 
nodes. This range is found to be common for all the sys¬ 
tems we have analysed for a fixed M. Thus, the primary 
criterion provides a small uniform range of e for many 
standard chaotic systems. 

In order to ensure that the RN is a true represen¬ 
tation of the chaotic attractor, we apply an additional 


FIG. 2: (Color online) Variation of the second smallest eigen 
value A 2 of the Laplacian matrix as a function of t for the 
RNs from Lorenz (filled circles connected by red solid line) 
and random (filled triangles connected by blue dashed line) 
time series. Results are shown for RNs with N = 2000 and 
5000, constructed with M = 3. The dotted line is a reference 
for zero. 


constraint that RN measures from standard chaotic time 
series are significantly different from that of a random 
time series. Though this condition appears to be sub¬ 
jective, we show below that it gives consistent results 
for all the standard chaotic systems considered in this 
work. For this, we compute the three important network 
measures, LD, CC and CPL for the RN from some stan¬ 
dard chaotic times series as a function of e and compare 
them with the corresponding measures from the RN of 
random time series. The equations for computing these 
measures are discussed in detail in the next section. We 
then find that the measure CPL is a good candidate to 
apply this additional constraint. This is shown in Fig. [3] 
for two standard chaotic time series. For e below 0.1, 
there are multiple disconnected networks with no giant 
cluster and the CPL computed is for the largest compo¬ 
nent. It is evident from the figure that there is a small 
range of e, say Ae, (marked by the two vertical dashed 
lines) where the difference in the value of CPL of the 
RN from chaotic time series and that of random is max¬ 
imum and as e increases above this range, the CPL for 
RN from all chaotic time series approaches that of ran¬ 
dom time series. More importantly, this range, which we 
call the critical range, is found to be approximately iden¬ 
tical for all the systems we have analysed and coincides 
with the common range found above corresponding to 
the emergence of the giant component in the RN. How¬ 
ever, it should be emphasized that this range, Ae, is an 
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FIG. 3: (Color online) Characteristic path length as a func¬ 
tion of threshold for RNs generated from Lorenz (filled circles 
connected by red solid line at the top), Rossler (filled squares 
connected by blue solid line in the middle) and random (filled 
triangles connected by dotted line) time series with N = 2000 
and M = 3. The two vertical dashed lines indicate the criti¬ 
cal range of e where CPL from chaotic time series differs from 
that of random time series. The critical threshold e c is the 
minimum of this range. Below e c , there is no giant cluster. 


empirical result and hence it is difficult to set any specific 
criterion for the upper bound either in terms of M or N. 
We choose the minimum value of this range as the critical 
threshold e c , for the construction of RN which we believe 
will capture the characteristic properties of the attractor. 
Nevertheless, we have checked and confirmed that any e 
within Ae does not make any qualitative change in the 
degree distribution of the RN and the related network 
measures to be presented below. Note that we do not 
follow the condition A 2 > 0 strictly as chosen by Eroglu 
et al. [35| , as this makes e c slightly different for different 
RNs. We have constructed the RN using the Gephi soft¬ 
ware (https://gephi.org/) taking time series from several 
standard chaotic attractors for e ranging from 0.05 to 
0.25, taking M = 3. We find that there is no giant clus¬ 
ter for e < 0.1 while the network becomes over connected 
for e > 0.15, in all cases. 

We now consider how the value of e c changes with M 
and N. To study this, we generate RNs from a number of 
low dimensional chaotic systems, both discrete and con¬ 
tinuous, by varying N from 2000 to 10000 and M from 
2 to 5. For each N and M, we scan a range of e val¬ 
ues between 0.02 to 0.20. The results from this detailed 
numerical analysis are compiled in Fig. [I] separately for 
continuous and discrete systems. The variation of < l > 
corresponding to Ae as N and M changes is shown with 


FIG. 4: (Color online) The figure shows how the critical range 
Ae (see text) varies with respect to N and M. The left panel 
shows Ae for two values of N as indicated, where the top 
panel is for continuous systems with M = 3 and the bottom 
panel for discrete systems with M — 2. The right panel shows 
the variation of Ae with M for N fixed at 10000. In the top 
panel, the red solid line is for Lorenz attractor and the blue 
dashed line is for Rossler attractor, while in the bottom panel, 
the same is for Lozi and Henon attractors respectively. In all 
figures, the dotted line is for random time series. 


that of random time series as reference. The left panel 
shows the dependence on N and the right panel shows the 
dependence on M. To study the dependence on N, we 
use the natural dimension of the system, namely, M = 3 
for continuous systems (top) and M — 2 for discrete sys¬ 
tems (bottom). Note that Ae for M = 2 has been shifted 
to 0.06 — 0.08. Moreover, in both cases, (M = 2 and 3), 
there is only a small decrease in Ae as N is increased 
from 2000 to 10000. This means that one can effectively 
use the same e for this whole range of N values. 

In the right panel, we show the effect of increasing the 
embedding dimension from the natural dimension of the 
system, with N fixed as 10000. Note that the value of e c 
clearly shifts with M which implies that each M requires 
a corresponding e c for the generation of the RN. But the 
more interesting result is that, for all the systems that 
we have analysed, we are able to find approximately the 
same critical range Ae corresponding to each M. This 
range is given in Table 1 where, e c in the third column 
is the critical threshold used by us corresponding to each 
M for all the further computations in this paper. Thus, 
even though our criterion for the selection of e c is not 
completely novel, we are able to provide a certain level of 
non-subjectivity in its choice which, we hope, will make 
the application of RN analysis more effective. 
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It may be noted that the above choice of the critical 
range of e is, in fact, analogous to and motivated by the 
selection of a scaling region in the conventional nonlinear 
time series analysis for deriving dynamical invariants like 
D 2 . We have already shown that the scaling region for 
many low dimensional chaotic systems can be selected 
algorithmically for a non-subjective computation of D 2 
and K 2 [3fl |39|. It is well known that the choice of the 
scaling region critically depends on the embedding di¬ 
mension M , while the change is not much if the number 
of data points changes from 2000 to 10000. In this way, 
the above result of getting an approximately same range 
of Ae for different chaotic systems is not very surprising. 


M 

Critical Range of e 

€c 

2 

0.06 - 0.08 

0.06 

3 

0.10 - 0.13 

0.10 

4 

0.14 - 0.18 

0.14 

5 

0.18 - 0.22 

0.18 


TABLE I: Critical range of e obtained empirically for each 
value of M 

We now give a simple mathematical explanation for the 
observed numerical results regarding the relation between 
e and M. Consider a random distribution of N points in 
M dimension. After uniform deviate transformation, the 
volume of the embedding space is unity and the average 
density of points < p >= N. The average separation 
between two points along any direction is 

< d >^(l) 1 / M^(_l_) 1/ M (4) 

This gives the critical value of distance for a given M 
below which the degree of a node tends to zero on the 
average. Hence e c for the random RN for each M must 
be sufficiently greater than < d >. 

To get more insight on the result given in Table 1 and 
know how e c varies for higher M values, we consider the 
limiting value of e, say e/, at which the RN is fully con¬ 
nected. That is, every node is connected to every other 
node with the degree of each node (N — 1) and the LD 
reaches its maximum possible value 1. Since we consider 
a uniform deviate, the size of the attractor is unity. For 
M = 2, the value of e at which this happens is the di¬ 
agonal length of the square, that is, e/(2) = \[2. For 
M = 3, e/(3) increases to >/3 and one can easily show 
that, in general, e/(M) = \/M. Note that this result 
is independent of N. This also implies that the differ¬ 
ence between successive e/, e/(M) — e/(M — 1), slowly 
decreases with M. Since the LD corresponding to e c is 
effectively a fraction of the total LD, one expects roughly 
the same dependence for e c on M, that is e c oc \[M. We 
have numerically verified this for the random RN for M 


upto 7. This means that, though the difference between 
successive e c appears to be almost a constant for small 
M as given in Table 1, this difference decreases slowly as 
M increases. 

As a final test to validate our choice of e c , we under¬ 
take a counter check by computing D 2 of some standard 
chaotic attractors using the RP (which is equivalent to 
the adjacency matrix) corresponding to e c . The method 
proposed by Thiel et al. [30| is used for this purpose. 
This method uses the cumulative probability distribu¬ 
tion p c (l) of the diagonal lines in the RP, corresponding 
to two different thresholds e and e+Ae using the relation: 



We have found that the D 2 values obtained in all cases 
are closer to the standard values for e = e c . This implies 
that the geometric complexity of the attractor is truly 
reflected in the RN constructed by our scheme. 


IV. MEASURES FROM RECURRENCE 
NETWORK 

We now compute the important network measures 
from the RN of several low dimensional chaotic attrac¬ 
tors, including discrete systems (maps) in two dimensions 
and continuous systems (flows) in three dimensions. For 
all systems, M is varied from 2 to 5 using the e c corre¬ 
sponding to each M with N varied from 2000 to 10000. 


A. Characteristic Path Length, Clustering 
Coefficient and Link Density 

We first compute the CPL, CC and LD from the RN 
and study their dependence on N and M. The equa¬ 
tions for computing these measures have been discussed 
in detail in the literature SI- The CPL is given by the 
equation 

1 N 

<l>= N(N^l)^ ls (6) 

where l s is the shortest path length for all pair of nodes 
{i,j) in the network. The maximum value of l s is taken 
as the diameter of the network, Id- If ki is the degree of 
the i th node, then 

1 N 

LD =m^7)^- ki (7) 

The CC of the network is defined through a local cluster¬ 
ing index c v . Its value is obtained by counting the actual 
number of edges in a sub graph with respect to node v 
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FIG. 5: (Color online)Variation of CPL (top panel) and CC 
(bottom panel) with N for RNs generated from Lorenz (filled 
circles connected dotted line), Rossler (filled triangles con¬ 
nected by solid line) and random (dashed line) time series. 
The value of M is fixed as 3. 


as reference to the maximum possible edges in the sub 
graph: 


Sij AiAjA 


IJ^JV 


k v {k v 1 ) 


( 8 ) 


The average value of c v is taken as the CC of the whole 
network: 


cc = ^J2 c » ( 9 ) 

V 

In Fig. [5] top panel, we show the variation of < l > with 
N for two standard chaotic systems for M = 3. RN from 
random time series is also added for comparison. Note 
that in all cases, < l > initially decreases with N, but 
saturates as N —> 10000, with the value of < l > for RN 
from random time series always less compared to that 
from chaotic systems. The variation of < l > with N 
can be understood from the degree distribution of the 
RN discussed in detail in the next section. We find that 
as N increases, the average degree of the nodes < k > 
also increases correspondingly. Typically, as N increases 
to 2IV, < k > shifts approximately to < 2k >, reducing 
< l >. In the bottom panel, we show the variation of CC 
with N which is found to be approximately constant for 
all systems for the range of N used. 

We have also studied numerically the variation of < l > 
and CC with M by fixing N = 10000 and the results 
for the same two systems are shown in Fig. [6j While 
CC remains constant, the variation of < l > is more 


FIG. 6: (Color online)Same as the previous figure, but the 
variation is with respect to M for a fixed N = 10000. Note 
that CPL saturates at the natural dimension of the system. 


interesting, showing saturation for M equal to or greater 
than the actual dimension of the system. We will show 
below that the same is true for degree distribution as well, 
leading to an important practical application of network 
measures. 


B. Degree Distribution 


We now consider the most important measure of a net¬ 
work, namely, the degree distribution P(fc) versus k. In 
Fig. [3 we show the degree distribution of the RN from 
Lorenz and Rossler attractor time series for N = 5000 
and 10000 by using M = 3 and e = 0.1. Note that the 
error bar is estimated from counting statistics resulting 
from the finiteness in the number of nodes. The statisti¬ 
cal error associated with any counting of n(k) is n(k). 
If n(k) —> 0, one typically takes the error to be normalised 
as 1. Thus, the error associated with P(k) is typically 

\/n(k) 

v ' and becomes 1/iV as n(k ) —I 0. It is evident from 
the figure that as N is doubled, the degree k of each node 
gets approximately doubled resulting in a shift along the 
X-axis and the range [k m i n ,k max \ of k values is shifted 
approximately to twice the range. Correspondingly, the 
P(k) values are reduced since the area under the distri¬ 
bution curve is constant. Assuming that the k values are 
continuous within the range [ k m i n , k ma x], we can write 

J2^kP(k) = l ( 10 ) 

k 























where A k(= 1) is the difference between successive k 
values. 

Changing k to k = jj, A k = ^ and P(k ) changes 
to P(k') so that 

^ Ak P(k) = 1 (11) 

k' 

Substituting for Afc in Eq.(10), we get 

P(k') = NP(k) = n(k) (12) 


It is convenient to represent the degree distribution in 
the rescaled variables as shown in Fig. [5] Note that 
the degree distribution for the two N values have now 
become identical and can be considered as almost sta¬ 
tionary in the rescaled variable k/N apart from small 
statistical fluctuations. To capture the real trend in the 
distribution, we show in Fig. [9] the same distributions 
shown in the previous figure without error bar and the 
values connected by line. 

One expects the scale invariance for the RN from a 
random time series whose degree distribution is Poisso- 
nian which can be approximated as Gaussian for large N. 
This is because, the degree of every node increases by an 
average value as N increases and the degree distribution 
appears identical in the rescaled variable k/N as can be 
seen from Fig. 1101 Here we find that an arbitrary de¬ 
gree distribution from the RN of a chaotic attractor also 
shows this property. A possible explanation, for attrac¬ 
tors whose measure is continuous is that, as the dynam¬ 
ical system evolves, the structure of the attractor also 
evolves in such a way that the probability density over 
the attractor is preserved once the basic structure of the 
attractor is formed. The degree ki of a reference node i 
represents the local connectivity of the RN and it corre¬ 
sponds to the local phase space density around the refer¬ 
ence point in the chaotic attractor from which the RN is 
constructed. It is well known that the local phase space 
density of the chaotic attractor is preserved in the RN 
[Tbt . Considering an infinitesimal hyper volume 14f(f) in 
M-dimension with radius e about a reference point f) in 
phase space, one can write jTif] : 


■4(fci(<0 + 1) ~ [ p(r)d M r^V M (e)p(fi) (13) 


where p(r)) is the invariant density around f). Note that 
in LHS, 1 is added to include the reference node (self 
loop). This gives a relation between the local measure in 
an attractor and that of a RN: 


p(ri ) = lim lim 

e—>0 N—too 


{ki + 1 ) 
V M {e)N 


(14) 


The above equation tells us that for the RN constructed 
with e c , the local probability density around a point on 
the attractor gets mapped to the degree of the corre¬ 
sponding node in the constructed RN. Since every point 


on the attractor is converted to a node on the RN, 
points with the same probability density will correspond 
to nodes with the same degree. The degree distribu¬ 
tion tells how many nodes have a given degree in the 
RN. This is equivalent to finding how many local regions 
on the entire attractor have the same probability den¬ 
sity. As we change from the phase space domain of the 
attractor to the domain of the network, the degree dis¬ 
tribution represents the global statistical measure of the 
probability density variations over the entire attractor. 
Thus, it seems natural that the degree distribution of 
the RN from any chaotic attractor shows the scale in¬ 
variance. The small deviations in the degree distribution 
as N increases is the result of the corresponding small 
fluctuations in the probability density. Also, the range 
of k values in the RN is a measure of the range of varia¬ 
tion of p(r) over the attractor. However, a direct relation 
connecting the probability distribution over the attractor 
and the degree distribution of the RN seems to be highly 
nontrivial owing to the fractal geometry of the attractor. 

From the above discussion, it becomes clear that a 
peak at high k value near k max in the degree distribution 
implies large number of relatively dense regions of high 
probability over the entire attractor. Many peaks in the 
degree distribution are indicative of large local density 
fluctuations over the attractor. For example, from Fig. [8] 
and Fig. [9l the fluctuation is much large for the RN from 
the Lorenz attractor compared to that of the Rossler at¬ 
tractor, though both have approximately the same range 
of [ kmin, kmax}- In this sense, one can say that the degree 
distribution is a characteristic measure of the structural 
complexity of an attractor. Note that this idea has been 
pointed out by many authors before [HI, |42| . Once the 
basic structure of the attractor is formed, a further in¬ 
crease in the number of nodes does not change the degree 
distribution qualitatively. In other words, RN analysis 
appears to be a useful tool to get meaningful results re¬ 
garding structural and topological properties of the at¬ 
tractors with less number of data points. We now show 
that there is a part in the degree distribution that cor¬ 
responds to the Poisson distribution where, the k values 
occur more by chance than by choice. 

For a random time series embedded in M-dimensional 
space, after being converted into uniform deviate, the av¬ 
erage density of points < p >= N. Hence the average 
number of points inside a M-dimensional sphere of ra¬ 
dius e is k ran = NVm , where Vm is the volume of the 
sphere. When the time series is converted into a RN, the 
condition for two nodes to be connected is that the dis¬ 
tance is < e. In other words, for random RN, typically a 
node is connected to k ran other nodes. Or, most nodes 
will have degree k ran and the degree distribution tends 
to be a Poissonian around this value. 

Now, for a non random time series, there will be sig¬ 
nificantly more nodes with degree greater than k ran and 
it is these nodes which describe the structure of the un¬ 
derlying attractor. Nodes with degree ~ k ran occur more 
by chance association rather than the true description of 
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FIG. 7: (Color online) Degree distribution of the RN gen¬ 
erated from the Lorenz attractor (top panel) for N = 5000 
(open red circles) and N = 10000 (filled green triangles ap¬ 
pearing in light gray shade in print), as indicated. Corre¬ 
sponding results for Rdssler attractor are shown in the bottom 
panel. In both cases, we use M = 3 and e c = 0.1. 



k/N 


FIG. 9: (Color online) The same distributions shown in the 
previous figure without error bar and connected by line, for 
clarity. The red solid line is for N = 5000 and the blue dashed 
line is for N = 10000. 

the system. Thus, characteristic information regarding 
the system is given by the nodes with degree > k ran and 
hence the value of k ran should be small compared to the 
range of k values in the degree distribution. 

Note that the position of k ran hr the degree distribu¬ 
tion depends on the choice of e and M. One expects a 
small peak around k ran in all degree distributions which 
becomes less significant as N increases. We now estimate 
the value of k ran for the choice of e and M that we use 
to compute the degree distributions of chaotic attractors. 
For M = 3, Vm = §7re 3 and hence 

kran ~ (15) 

With e = 0.1, we get p(p ~ 0.004, which is independent 
of N. This is sufficiently small compared to the rescaled 
jj values as can be seen from Fig. [5] and Fig. [9l where 
pp « 0.04. Since pp oc e 3 for M = 3, a small increase 
in e can shift the Poisson range significantly to the right 
in the degree distribution. This shows the importance of 
choosing the correct e c . 

For the discrete systems with M = 2, we have 


FIG. 8: (Color online) Rescaled distributions (see text) of 
the Lorenz (top panel) and Rdssler (bottom panel) attractors 
shown in the previous figure for N = 5000 (open red circles) 
and N = 10000 (filled green triangles appearing in light gray 
shade in print). Note that, after rescaling, the distributions 
for the two N values become statistically identical in both 
cases. 


N 


7T£ 


(16) 


Using the optimum value of e = 0.06 used for M = 2, we 
have ss 0.011 which is << ppS as will be shown 
below for discrete systems. Finally, for M = 4, we have a 
hyper cube of unit volume. The general formula for the 
volume of a Euclidean ball of radius e in M-dimension 
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FIG. 10: (Color online) Degree distribution of the RN from 
the Rossler attractor time series (filled blue triangles) in com¬ 
parison with that of the random time series (filled red circles 
seen as Poisson distribution on the left). The top panel is for 
N = 5000 and the bottom panel for N = 10000, with M = 3 
and e = 0.1 in both cases. We show the rescaled distribu¬ 
tions so that both panels appear identical. The two vertical 
lines indicate the Poisson statistics part of the distribution 
or Poisson statistic error (denoted PSE) which coincides with 
the degree distribution of the random time series. It is shown 
magnified in the inset. 


(for even M) is: 


V M (e) 


7T M / 2 

(¥)! 


M 


(17) 


For M = 4, V±{e) = For the threshold e = 0.14, we 
find ^ « 0.00196. 

The above results are explicitly shown in Fig. [TO] and 
Fig. |TT] In Fig. [TUI we show the rescaled degree distribu¬ 
tions of RNs from random and Rossler attractor time se¬ 
ries plotted together for N = 5000 and 10000. Note that 
the Poisson distribution part, shown by the two vertical 
lines, almost exactly coincides with the degree distribu¬ 
tion of the random time series. This part is shown magni¬ 
fied in the inset. In Fig. [TlJ we show the rescaled degree 
distributions of the Cat map and the random time series 
together for M = 2 (top panel) and for M = 3 (bottom 
panel). For M = 2, both the distributions are almost 
identical and peak exactly at k ra n = 0.011 in agreement 
with our calculations above. For M = 3, the peak for 
the random distribution is shifted to 0.004 as expected, 
while that for Cat map remains almost unchanged and 
hence both the distributions can be easily differentiated. 

We have, so far, computed the degree distribution by 
taking the actual dimension of the attractor. However, in 


FIG. 11: (Color online) Degree distributions of the RNs gen¬ 
erated from Cat Map (filled blue circles) and random (filled 
red triangles) time series for N = 5000 for two values of M, 
as indicated. Note that, in going from M = 2 to M = 3, the 
distribution of the random network is shifted while that for 
the Cat Map remains stationary. 
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FIG. 12: (Color online) Rescaled degree distributions (with¬ 
out error bar for clarity) for the RNs from Lorenz (top 
panel) and Rossler (bottom panel) attractor time series for 
M = 3(open red circles), M = 4(Hlled blue triangles) and 
M = 5(black star like points). 
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FIG. 13: (Color online) Degree distributions (without error 
bar) of the RN constructed from the time series of Chen 
hyperchaotic attractor for M = 2 (thick blue dotted line), 
M = 3 (thick green solid line), M = 4 (black dashed line) 
and M = 5 (thin red solid line). For each M, the correspond¬ 
ing value of e c given in Table 1 is used for constructing the 
RN with N = 5000. 


the analysis of the real world data, there is no a priori in¬ 
formation regarding the dimension of the system. Then 
one has to compute the network measures for different 
M values and check for saturation. In Fig. [T3] we show 
the rescaled degree distributions of RNs from Lorenz and 
Rossler attractors for M = 3,4, 5. For each M, the cor¬ 
responding e c value found empirically, as given in Table 
1, is used. We find that the degree distribution remains 
almost invariant (apart from small changes due to the ef¬ 
fect of embedding), for M values equal to or greater than 
the actual dimension of the system. On the one hand, it 
is a counter check whether we are using the correct value 
of e c for each M and on the other hand, the result tells us 
that the usual practice of using a high value of M for real 
world data works for network measures as well provided 
the correct e c corresponding to each M is used. How¬ 
ever, a higher M requires a correspondingly larger value 
of N. We find that it is sufficient to use N < 10000 and 
M < 5 for a proper characterization of low dimensional 
chaotic systems using RN measures. This also leads to 
a practical application of the proposed method discussed 
in more detail in the next section. 

We next show that the present scheme is efficient not 
only for low dimensional chaotic attractors, but also for 
the analysis of high dimensional and hyperchaotic sys¬ 
tems as well. In Fig. [13l we show the degree distribu¬ 
tions of the RN constructed from the time series of Chen 
hyperchaotic flow [43| for M varying from 2 to 5, with 


corresponding values of e c . We generate the hyperchaotic 
time series with 5000 data points using the standard set 
of parameters: a = 35, b = 4.9, c = 25, d = 5, e = 35 and 
k = 100. Note that the degree distributions for M = 4 
and 5 are almost identical while that for lower M values 
deviate, since the attractor dimension is > 3. This is also 
a direct confirmation of our argument above that for a 
real world data whose dimension is unknown, one has to 
check for saturation of network measures by increasing 
the M value. This is somewhat equivalent to finding the 
saturated Di value by increasing M in the conventional 
nonlinear time series analysis. 

Another important outcome of our scheme is that we 
are able to compare the characteristic measures derived 
from the RN of different chaotic attractors since the 
analysis is done using identical threshold for a fixed M. 
For example, the degree distribution of the RN typi¬ 
cally characterises the structure of the attractor, as dis¬ 
cussed above. Hence, a visual inspection of the degree 
distribution can provide some qualitative information on 
the structural complexities of standard low dimensional 
chaotic attractors, as shown in Fig. [14] and Fig. [I5l The 
degree distribution in each case is the average from four 
RNs generated using different initial conditions for the 
attractor. From a comparison of the degree distributions 
in Fig. [14] one can infer that the Lorenz attractor is struc¬ 
turally more complex compared to the other three since 
the fluctuation in the probability density is maximum for 
it. More accurate results can be obtained from a quanti¬ 
tative analysis of the network measures. 

For the logistic map, we use the fully chaotic region 
with M = 1 and e c = 0.01 and hence the very large peak 
in the degree distribution corresponding to that e value 
is due to Poisson statistics. The logistic map requires a 
special mention. For an attractor in one dimension, the 
Poisson value of k/N = e c , the threshold itself. In other 
words, for a random distribution in the unit interval, the 
degree distribution is typically a Poissonian around de¬ 
gree k/N = e c in our scheme. However, for the RN from 
the logistic attractor, depending on the probability den¬ 
sity variations, the number of degrees of any node can be 
>> e c or << e c . For example, from the figure, there are 
nodes with degree as high as 0.09 and as low as 0.002. 
The structure of the logistic attractor is actually charac¬ 
terized by jt values > e c . 


V. APPLICATIONS 

So far, we have been discussing the construction and 
analysis of RN from chaotic time series. It is also impor¬ 
tant to know how effectively our scheme can be applied 
to time series data from the real world. An important 
difference is that for standard chaotic systems, the di¬ 
mensionality of the system is known a priori and Ad can 
be fixed accordingly while for real world data, this in¬ 
formation is absent. In the conventional nonlinear time 
series analysis, one computes dynamical invariants as a 
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FIG. 14: (Color online) Characteristic degree distributions of 
four standard low dimensional chaotic attractors for M = 3. 
The systems are (a) Lorenz (b) Rossler (c) Duffing and (d) 
LTeda attractor. Standard parameter values given in ficj are 
used for the generation of time series from Duffing and Ueda 
attractors. The average degree distribution for RNs generated 
from four initial conditions is shown in all cases. 



FIG. 15: (Color online) Same as the previous figure for four 
standard chaotic maps, namely, (a) Henon (b) Lozi (c) Cat 
Map and (d) Logistic Map. For the logistic map in the fully 
chaotic regime, we use M = 1 and e c = 0.01 which results 
in the high peak corresponding to that e value. For all other 
cases, we use M = 2. 


function of M and check for saturation with respect to 
M. For RN analysis, what is normally done is to use a 
sufficiently large value of M to ensure a proper embed¬ 
ding of the underlying attractor jlj)j. Our scheme indi¬ 
cates that a large M requires a correspondingly large N 
and e c . However, the number of data points in real world 
time series is normally less (say, < 10000) and contami¬ 
nated by different types of noise. One of the advantages 
of RN analysis is that it is possible to get information 
regarding the underlying system from a time series with 
limited number of data points. Hence for practical imple¬ 
mentation of the scheme for real world data, we suggest 
that it is better to start the computation taking a small 
M, go for higher dimensions successively and check for 
saturation of the network measures for two successive 
dimensions, rather than using a very high embedding di¬ 
mension to start with. 

To illustrate the potential applications of our ap¬ 
proach, we consider two examples. In the first example, 
we show that the proposed scheme is capable of iden¬ 
tifying the dimensionality of the underlying system and 
the presence of white noise in real world data. For this, 
we present the RN analysis of the light curves from a 
black hole system GRS 1915+105. The light curves from 
this black hole system have been classified into 12 spec¬ 
troscopic states by Belloni et al. [4l| and we take light 
curves from two representative states 9 and y, which are 
shown in Fig. 1161 In an earlier paper [37t | . we have shown 
by computing D 2 that the state 9 has signatures of de¬ 
terministic nonlinear behavior (with D 2 < 3) and y is 
white noise. We construct the RN from the two time se¬ 
ries for different values of M starting from M = 2 using 
the e c corresponding to each M and compute the network 
measures in each case. We find that if the system is of 
finite dimension, the measures saturate beyond a certain 
M which is taken as the dimension of the system. In 
Fig. E3 we show the rescaled degree distributions of the 
two light curves for M = 2,3 and 4. Note that the degree 
distributions for the two states are completely different. 
For the y state, they are Poissonian and shows the typical 
shift without any saturation as M increases, indicating 
pure white noise. On the other hand, the state 9 is qual¬ 
itatively different with the degree distributions getting 
saturated for M = 3 and 4 and remains stationary for 
any higher embedding dimension. 

Due to the inherent non-subjectivity in the choice of e c , 
the scheme is also ideal for the surrogate analysis using 
network measures CC and CPL as discriminating statis¬ 
tic to detect deterministic nonlinearity in real world data. 
It can be used as complementary to conventional anal¬ 
ysis with measures like D 2 and K 2 and has the added 
advantage that the length of the data required is much 
less compared to conventional methods. 

As the second application, we consider detection of a 
dynamical transition using the RN measures derived from 
our method. This is known to be an important applica¬ 
tion of RNs pH EH- The example we choose is the chaos- 
hyperchaos transition in a time delayed system, namely, 
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FIG. 16: (Color online) Part of the astrophysi- 

cal light curves from the black hole system GRS 
1915 + 105 for two temporal states, 9 and \ (Ref: 
xte.mit.edu/ ehm/1915_lightcurves.html). 


FIG. 18: (Color online) Variation of the GPL and CC of the 
RN constructed by our scheme from the time series of M-G 
system for the range of td values representing the transition 
from chaos to hyperchaos. Corresponding values of the second 
largest LE (denoted LE 2 ) are shown in the top panel. The 
dashed vertical line indicates the transition point. 



k/N 


FIG. 17: (Color online) Rescaled degree distributions com¬ 
puted from the RNs for the two light curves 9 and \ for 
three M values. The distributions are for M = 2 (filled black 
squares), M = 3 (open red circles) and M = 4 (filled green 
triangles appearing in light gray shade in print) for 9 state. 
The distributions for the three M values for the \ state can be 
clearly distinguished, as indicated. In both cases, N = 3000. 
Note that for the 9 state the two degree distributions for 
M = 3 and 4 almost coincide. 


the Mackey-Glass (M-G) system |45j given by the equa¬ 
tion: 


dx /3 x td 
dt 1 + x™ D 


(18) 


The constants /3, 7 , To and n are real numbers and x TD 
represents the value of the variable x at time t — to- 
Depending on the values of the parameters, the system 
displays a range of periodic, chaotic and hyperchaotic 
dynamics. We fix the parameters /3 = 2 , 7 =l,n =10 
with To as the control parameter. As the value of the 
parameter to increases, the asymptotic state of the sys¬ 
tem changes from periodic to chaotic and then to hyper¬ 
chaos. We have recently shown [4f| using a dimensional 
analysis that the transition to hyperchaos occurs at a 
critical value to = 4.038. Here we check whether the 
RN measures obtained from our scheme can detect this 
transition. 

To do the analysis, we first generate time series of 
length 100000 from the system for to values ranging 
from 3.5 to 4.5 increasing in steps of 0.1. Each of this 
time series is then split up into 10 time series of length 
Nt = 10000. We construct RN from all these and com¬ 
pute the measures CC and CPL taking M = 4 and 
e = 0.14, which is the minimum M value required for 
a hyperchaotic attractor (any higher embedding dimen¬ 
sion is equally good). Using the TISEAN package j47j . 
we also compute the second largest Lyapunov Exponent 
(LE) (denoted LE 2 ) which crosses zero as the system 
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cross over to liyperchaotic phase. The results of compu¬ 
tations are shown in Fig. [THl where the vertical dashed 
line indicates the transition point. The error bar comes 
from the standard deviation of the values obtained from 
the ten time series used for each td- It is clear that the 
present scheme can detect the transition as both CC and 
CPL show a discontinuity at the transition point. 

VI. DISCUSSION AND CONCLUSION 

Recurrence networks have become important statis¬ 
tical tools for characterising the structural properties 
of chaotic attractors. They are complex networks con¬ 
structed from chaotic time series using a suitable scheme 
that maps information inherent in the time series into the 
network domain. One can then use the network based 
measures to quantify the geometric properties of the un¬ 
derlying attractor. The distinct advantages of such an 
exercise are that the network measures can be efficiently 
computed from less number of data points in the time se¬ 
ries and these measures can also be used as complimen¬ 
tary to conventional measures of nonlinear time series 
analysis. Here we present a general scheme to construct 
the RN from a chaotic time series that can be applied 
equally well to time series from standard chaotic attrac¬ 
tors as well as real world data. We use identical value of 
threshold to construct the RN from different time series 
for a given embedding dimension. The scheme is thus 
found to be suitable to compare the structural proper¬ 
ties of two chaotic attractors by computing the statistical 
measures of the corresponding recurrence networks. To 
illustrate how the scheme can be implemented in prac¬ 
tice, it is used for the analysis of light curves from a black 
hole system and to identify the transitions between two 
dynamical regimes in a time delayed system. 

It is important to know the merits, limitations and po¬ 
tential applications of the proposed scheme for its proper 
implementation. Several methods have already been pro¬ 
posed in the literature for the construction of RN, as dis¬ 
cussed in §1. These methods mainly differ in the criteria 
for the selection of e c . The main aspects in our approach 
compared to the earlier methods are the uniform devi¬ 
ate transformation of the time series and the criterion 
used for the selection of e c with its linkage to embedding 
dimension M. These changes make it possible to look 
for a uniform critical threshold for different chaotic time 
series. 

We do not claim that the RN constructed by our ap¬ 
proach is optimum for all the different types of time series 
and applications. The critical range Ae presented in Ta¬ 
ble 1 for each M is an empirical choice resulting from the 
analysis on a limited number of nodes (N < 10000). It 
is primarily motivated to get a uniform range for differ¬ 


ent chaotic attractors as explained in §111 and is not a 
rigorous result. Though we expect the range to be valid 
in general, it may require improvement for specific time 
series and more accurate applications. However, by get¬ 
ting an identical value of e c for different time series, we 
are able to achieve a certain level of non subjectiveness in 
the construction of RN, especially for the analysis of time 
series from the real world, as we have shown explicitly. 

Another important outcome of the present approach 
not reported previously, is the realisation that the value 
of e c should be linked to M. This implies that the choice 
of M is equally important for RN construction from 
time series and a very large value of M, as is generally 
believed, may not provide optimum result with limited 
number of data points. It is also interesting to note that 
there is always a part in the degree distribution of the RN 
from a chaotic attractor that corresponds to Poisson dis¬ 
tribution where the k values occur more by chance than 
by choice. Its position in the degree distribution depends 
on the choice of e c and M. 

There are at least four important potential applica¬ 
tions for our approach, three of which we have shown 
here explicitly: 

i) to compare the structural properties of two chaotic 
attractors using network measures through the construc¬ 
tion of RN 

ii) to study the transition between two dynamical 
regimes as a control parameter is varied 

iii) to identify the dimension of the underlying attrac¬ 
tor from the analysis of time series data and 

iv) for surrogate analysis using any of the RN mea¬ 
sures as discriminating statistic where a non subjective 
comparison of the network measures from data and sur¬ 
rogates is required. 

Finally, an important step forward in our analysis is 
to try and develop a similar scheme for RNs where the 
connections have weight factors. Here we have considered 
unweighted RNs so that the resulting adjacency matrix is 
binary. We hope that a weighted RN can unravel more in¬ 
formation regarding the topological and structural prop¬ 
erties of chaotic attractors. Another possible application 
of the scheme, that is important in the analysis of real 
world data, is to study the effect of noise on RN and the 
measures derived from it. These works are currently in 
progress and will be reported elsewhere. 
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